Modelo

fitSM1_Mensual <- tslm(tsSM1_Mensual ~ trend + season) 
summary(fitSM1_Mensual)

Call:
tslm(formula = tsSM1_Mensual ~ trend + season)

Residuals:
     Min       1Q   Median       3Q      Max 
-20522.7  -5190.9    215.7   6726.8  16702.8 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  66531.9     7294.5   9.121 4.21e-09 ***
trend         -173.1      200.7  -0.862 0.397393    
season2     -19137.6     9635.6  -1.986 0.059061 .  
season3      -1263.5     9641.8  -0.131 0.896879    
season4     -15856.4     9652.3  -1.643 0.114034    
season5      -6134.4     9666.9  -0.635 0.531966    
season6      -9988.0     9685.6  -1.031 0.313160    
season7     -27046.6     9708.4  -2.786 0.010505 *  
season8     -38674.8     9735.4  -3.973 0.000602 ***
season9     -10641.7     9766.3  -1.090 0.287159    
season10    -15040.7     9801.3  -1.535 0.138536    
season11     -6105.3     9840.3  -0.620 0.541071    
season12     -3378.2     9883.2  -0.342 0.735595    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 11800 on 23 degrees of freedom
Multiple R-squared:  0.5796,    Adjusted R-squared:  0.3603 
F-statistic: 2.643 on 12 and 23 DF,  p-value: 0.02184

Resultados:

• Valor del coeficiente de determinación: 0.5796. Es un valor bajo, el modelo únicamente explica el 58% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 20 días

forecastfitSM1_Mensual <- forecast(fitSM1_Mensual, h=12, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM1_Mensual, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Modelo

fitSM2_Mensual <- tslm(tsSM2_Mensual ~ trend + season) 
summary(fitSM2_Mensual)

Call:
tslm(formula = tsSM2_Mensual ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-26325  -4033   1223   6470  16873 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  83389.3     7617.6  10.947 1.35e-10 ***
trend         -898.9      209.6  -4.289 0.000274 ***
season2     -13475.5    10062.5  -1.339 0.193598    
season3       7508.4    10069.0   0.746 0.463406    
season4     -15023.1    10079.9  -1.490 0.149708    
season5      -7727.5    10095.2  -0.765 0.451777    
season6     -11683.0    10114.7  -1.155 0.259931    
season7     -19075.8    10138.6  -1.882 0.072621 .  
season8     -28408.6    10166.7  -2.794 0.010305 *  
season9      -9204.4    10199.1  -0.902 0.376160    
season10      4809.1    10235.6   0.470 0.642890    
season11     -1696.0    10276.3  -0.165 0.870356    
season12     -3191.5    10321.1  -0.309 0.759941    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 12320 on 23 degrees of freedom
Multiple R-squared:  0.6541,    Adjusted R-squared:  0.4736 
F-statistic: 3.624 on 12 and 23 DF,  p-value: 0.003891

Resultados:

• Valor del coeficiente de determinación: 0.6541. Es un valor bajo, el modelo únicamente explica el 65% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM2_Mensual <- forecast(fitSM2_Mensual, h=12, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM2_Mensual, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Modelo

fitSM3_Mensual <- tslm(tsSM3_Mensual ~ trend + season) 
summary(fitSM3_Mensual)

Call:
tslm(formula = tsSM3_Mensual ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-87112 -18443   4781  17647  82107 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   299382      23189  12.911 5.06e-12 ***
trend           1866        638   2.926  0.00761 ** 
season2       -50805      30631  -1.659  0.11077    
season3       -29034      30651  -0.947  0.35336    
season4       -64190      30684  -2.092  0.04768 *  
season5       -54165      30730  -1.763  0.09125 .  
season6       -85197      30790  -2.767  0.01097 *  
season7      -144850      30863  -4.693 1.00e-04 ***
season8      -171144      30948  -5.530 1.27e-05 ***
season9       -69624      31047  -2.243  0.03486 *  
season10      -52696      31158  -1.691  0.10430    
season11      -37281      31282  -1.192  0.24550    
season12        7092      31418   0.226  0.82341    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 37510 on 23 degrees of freedom
Multiple R-squared:  0.7577,    Adjusted R-squared:  0.6312 
F-statistic: 5.992 on 12 and 23 DF,  p-value: 0.0001238

Resultados:

• Valor del coeficiente de determinación: 0.7577. Es un valor que podría ser aceptable, el modelo únicamente explica el 76% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM3_Mensual <- forecast(fitSM3_Mensual, h=12, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM3_Mensual, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Submetering 1

### Descomponer la serie
#Componentes_SM1_Mensual <- decompose(tsSM1_Mensual)
### Plot
#plot(Componentes_SM1_Mensual )
## stl(tsSM1_Mensual)
### Check summary statistics for decomposed sub-meter 1
#summary(Componentes_SM1_Mensual )

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

 tsSM1_Mensual %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal              trend            remainder        
 Min.   :-17333.003   Min.   :46768.55   Min.   :-37011.52  
 1st Qu.: -5016.432   1st Qu.:48940.13   1st Qu.: -3448.01  
 Median : -3280.804   Median :51834.57   Median :  -575.02  
 Mean   :     0.000   Mean   :51611.28   Mean   : -1053.45  
 3rd Qu.:  5593.844   3rd Qu.:53148.20   3rd Qu.:  3802.57  
 Max.   : 14381.497   Max.   :57632.43   Max.   : 32860.67  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     10610         4208      7251         16690
   %  63.6         25.2      43.4         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.9001  0.9452  0.8181  0.9965  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 361 19 13
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 37 2 2
 $ inner: int 1
 $ outer: int 15

tsSM1_Mensual %>%
  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
  autoplot()

Remainder no es un muelle. Observamos tendencia creciente y fuerte componente estacional

Componentes_SM1_MensualSTL<- tsSM1_Mensual %>%stl( t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM1_MensualSTL, ts.colour = 'blue')

plot.ts(tsSM1_Mensual, col = 'gray')
lines(Componentes_SM1_MensualSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 2

## Descomponer la serie
Componentes_SM2_Mensual <- decompose(tsSM2_Mensual)
## Plot
plot(Componentes_SM2_Mensual )

## Check summary statistics for decomposed sub-meter 2
summary(Componentes_SM2_Mensual )

         Length Class  Mode     
x        36     ts     numeric  
seasonal 36     ts     numeric  
trend    36     ts     numeric  
random   36     ts     numeric  
figure   12     -none- numeric  
type      1     -none- character

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

 tsSM2_Mensual %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE) %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal              trend            remainder        
 Min.   :-21630.641   Min.   :48453.75   Min.   :-43864.84  
 1st Qu.: -4761.576   1st Qu.:48533.62   1st Qu.: -5720.49  
 Median :  1499.661   Median :56646.72   Median :  -144.02  
 Mean   :    -0.001   Mean   :59759.83   Mean   : -1096.91  
 3rd Qu.:  5474.915   3rd Qu.:70783.48   3rd Qu.:  3351.28  
 Max.   : 15664.732   Max.   :77538.14   Max.   : 20892.53  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     10236        22250      9072         23422
   %  43.7         95.0      38.7         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8448  0.9452  0.8851  0.9989  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 361 19 13
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 37 2 2
 $ inner: int 1
 $ outer: int 15

tsSM2_Mensual %>%
  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
  autoplot()

Componentes_SM2_MensualSTL<- tsSM2_Mensual %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM2_MensualSTL, ts.colour = 'blue')

plot.ts(tsSM2_Mensual, col = 'gray')
lines(Componentes_SM2_MensualSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 3

## ## Descomponer la serie
## Componentes_SM3_Mensual <- decompose(tsSM3_Mensual)
## ## Plot
## plot(Componentes_SM3_Mensual )
## ## Check summary statistics for decomposed sub-meter 2
## summary(Componentes_SM3_Mensual )

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia creciente

tsSM3_Mensual %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal              trend            remainder        
 Min.   :-107431.70   Min.   :249102.3   Min.   :-91715.37  
 1st Qu.: -12946.30   1st Qu.:250836.3   1st Qu.:-17758.23  
 Median :   9126.59   Median :276150.6   Median : -1364.69  
 Mean   :      0.00   Mean   :274156.0   Mean   : -2901.04  
 3rd Qu.:  25952.58   3rd Qu.:283886.9   3rd Qu.: 11936.12  
 Max.   :  68450.32   Max.   :316036.5   Max.   : 83277.64  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     38899        33051     29694         82030
   %  47.4         40.3      36.2         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8028  0.9452  0.8268  0.9889  0.9999 

 Other components: List of 5
 $ win  : Named num [1:3] 361 19 13
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 37 2 2
 $ inner: int 1
 $ outer: int 15

Componentes_SM3_MensualSTL<- tsSM3_Mensual %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM3_MensualSTL, ts.colour = 'blue')

plot.ts(tsSM3_Mensual, col = 'gray')
lines(Componentes_SM3_MensualSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Evolución mensual del consumo (sin estacionalidad). Años 2007-2009

#tsSM1_Mensual_Adjusted <- tsSM1_Mensual - Componentes_SM1_Mensual$seasonal
#autoplot(tsSM1_Mensual_Adjusted)

tsSM1_Mensual_detrended <- tsSM1_Mensual - Componentes_SM1_MensualSTL$time.series[,2]
plot.ts(tsSM1_Mensual_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM1_MensualSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM1_MensualSTL$time.series[,3])
qqline(Componentes_SM1_MensualSTL$time.series[,3])

shapiro.test(Componentes_SM1_MensualSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM1_MensualSTL$time.series[, 3]
W = 0.75719, p-value = 2.554e-06

## Volvemos a descomponer la serie
plot(stl(tsSM1_Mensual_detrended ,s.window = "periodic",robust = TRUE))

tsSM1_HW_Mensual <- HoltWinters(tsSM1_Mensual_detrended ,
## Holt Winters Exponential Smoothing & Plot
                              #  tsSM1_Mensual_Adjusted, 
  beta=TRUE, gamma=TRUE)

# optimiza
plot(tsSM1_HW_Mensual,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Mensual_for <- forecast(tsSM1_HW_Mensual, h=12)
plot(tsSM1_HW_Mensual_for , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Mensual_forC <- forecast(tsSM1_HW_Mensual, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Mensual_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM2_Mensual_Adjusted <- tsSM2_Mensual - Componentes_SM2_Mensual$seasonal
# autoplot(tsSM2_Mensual_Adjusted)

tsSM2_Mensual_detrended <- tsSM2_Mensual - Componentes_SM2_MensualSTL$time.series[,2]
plot.ts(tsSM2_Mensual_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM2_MensualSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM2_MensualSTL$time.series[,3])
qqline(Componentes_SM2_MensualSTL$time.series[,3])

shapiro.test(Componentes_SM2_MensualSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM2_MensualSTL$time.series[, 3]
W = 0.85783, p-value = 0.0002898

## Volvemos a descomponer la serie
plot(stl(tsSM2_Mensual_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_Mensual <- HoltWinters(tsSM2_Mensual_detrended ,
                                # tsSM2_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_Mensual,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_Mensual_for <- forecast(tsSM2_HW_Mensual, h=12)
plot(tsSM2_HW_Mensual_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Mensual_forC <- forecast(tsSM2_HW_Mensual, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Mensual_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

# ## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM3_Mensual_Adjusted <- tsSM3_Mensual - Componentes_SM3_Mensual$seasonal
# autoplot(tsSM3_Mensual_Adjusted)

tsSM3_Mensual_detrended <- tsSM3_Mensual - Componentes_SM3_MensualSTL$time.series[,2]
plot.ts(tsSM3_Mensual_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM3_MensualSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM3_MensualSTL$time.series[,3])
qqline(Componentes_SM3_MensualSTL$time.series[,3])

shapiro.test(Componentes_SM3_MensualSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM3_MensualSTL$time.series[, 3]
W = 0.93643, p-value = 0.03938

## Volvemos a descomponer la serie
plot(stl(tsSM3_Mensual_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_Mensual <- HoltWinters(tsSM3_Mensual_detrended ,
                              # tsSM3_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_Mensual,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Mensual_for <- forecast(tsSM3_HW_Mensual, h=12)
plot(tsSM3_HW_Mensual_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Mensual_forC <- forecast(tsSM3_HW_Mensual, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Mensual_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

RMSE_SM1_GranMensualFor<-accuracy(tsSM1_HW_Mensual_forC   ,DatosRealesSeriesMensuales$Sub_metering_1)[4]
RMSE_SM2_GranMensualFor<-accuracy(tsSM2_HW_Mensual_forC   ,DatosRealesSeriesMensuales$Sub_metering_2)[4]
RMSE_SM3_GranMensualFor<-accuracy(tsSM3_HW_Mensual_forC   ,DatosRealesSeriesMensuales$Sub_metering_3)[4]


ResEvolMensual2_SinEstacionalidadConTendencia <-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranMensualFor,
               RMSE_SM2_GranMensualFor,
               RMSE_SM3_GranMensualFor)
)
ResEvolMensual2_SinEstacionalidadConTendencia %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Evolución mensual del consumo (con estacionalidad). Años 2007-2009

Vamos a ver una predicción del consumo enegético de cada submedidor con el método de HW. En este caso, no eliminaremos la estacionalidad, para ver así si se mantiene la tendencia de consumo.

tsSM1_Mensual_AdjustedII <- tsSM1_Mensual # - Componentes_SM1_Mensual$seasonal
autoplot(tsSM1_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM1_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_MensualII <- HoltWinters(tsSM1_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Mensual_forII <- forecast(tsSM1_HW_MensualII, h=12)
plot(tsSM1_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Mensual_forCII <- forecast(tsSM1_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM2_Mensual_AdjustedII <- tsSM2_Mensual ## - Componentes_SM2_Mensual$seasonal
autoplot(tsSM2_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM2_Mensual_AdjustedII)) 

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_MensualII <- HoltWinters(tsSM2_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_Mensual_forII <- forecast(tsSM2_HW_MensualII, h=12)
plot(tsSM2_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Mensual_forCII <- forecast(tsSM2_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM3_Mensual_AdjustedII <- tsSM3_Mensual # - Componentes_SM3_Mensual$seasonal
autoplot(tsSM3_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM3_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_MensualII <- HoltWinters(tsSM3_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Mensual_forII <- forecast(tsSM3_HW_MensualII, h=12)
plot(tsSM3_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Mensual_forCII <- forecast(tsSM3_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

RMSE_SM1_GranMensualForII<-accuracy(tsSM1_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_1)[4]
RMSE_SM2_GranMensualForII<-accuracy(tsSM2_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_2)[4]
RMSE_SM3_GranMensualForII<-accuracy(tsSM3_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_3)[4]


ResEvolMensual2SinTendenciaConEstacionalidad<-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranMensualForII,
               RMSE_SM2_GranMensualForII,
               RMSE_SM3_GranMensualForII)
   
    
)

ResEvolMensual2SinTendenciaConEstacionalidad %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Comparación de resultados. Predicción antes y despues de eliminar estacionalidad

Comparacion<-cbind.data.frame(
  ResEvolMensual2_SinEstacionalidadConTendencia,
  ResEvolMensual2SinTendenciaConEstacionalidad[,2])

colnames(Comparacion) <- c("Submedidor","RMSE Con Estacionalidad","RMSE Sin estacionalidad")
Comparacion %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Resultados más fiables y con menor error trás haber eliminado la estacionalidad de la serie. También así podemos ver la tendencia.

#p1<-autoplot(tsSM1_Mensual_Adjusted)
#p2<-autoplot(tsSM1_Mensual_AdjustedII)
#grid.arrange(
#             grobs = list(p1,p2),
#           nrow = 2,
#             top="Sin estacionalidad",bottom="Con estacionalidad"
#             )

Modelo

fitSM1_Semanal <- tslm(tsSM1_Semanal ~ trend + season) 
summary(fitSM1_Semanal)

Call:
tslm(formula = tsSM1_Semanal ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-11708  -2978   1845   4714   7123 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  59656.6     4512.2  13.221 6.48e-09 ***
trend         -508.6      263.8  -1.928  0.07593 .  
season2      12293.0     5646.7   2.177  0.04850 *  
season3      23297.9     5665.1   4.113  0.00122 ** 
season4       7005.5     5695.7   1.230  0.24051    
season5      13062.8     5738.3   2.276  0.04039 *  
season6      81428.1     5792.6  14.057 3.07e-09 ***
season7      91174.1     5858.3  15.563 8.75e-10 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 6908 on 13 degrees of freedom
Multiple R-squared:  0.975, Adjusted R-squared:  0.9616 
F-statistic: 72.49 on 7 and 13 DF,  p-value: 2.122e-09

Resultados:

• Valor del coeficiente de determinación: 0.975. Es un valor bajo, el modelo únicamente explica el 97.5% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 20 días

forecastfitSM1_Semanal <- forecast(fitSM1_Semanal, h=7, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM1_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético semanal total del año 2010 \n  a través de un modelo de regresión")

Modelo

fitSM2_Semanal <- tslm(tsSM2_Semanal ~ trend + season) 
summary(fitSM2_Semanal)

Call:
tslm(formula = tsSM2_Semanal ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-34037 -13452  -3470   8762  38539 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  88194.9    15111.4   5.836 5.82e-05 ***
trend        -2641.6      883.3  -2.991 0.010427 *  
season2      43786.9    18910.7   2.315 0.037565 *  
season3      73582.1    18972.5   3.878 0.001903 ** 
season4      -4120.0    19075.0  -0.216 0.832350    
season5      28106.0    19217.7   1.463 0.167355    
season6      64649.2    19399.5   3.333 0.005399 ** 
season7      83987.4    19619.5   4.281 0.000895 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 23140 on 13 degrees of freedom
Multiple R-squared:  0.7717,    Adjusted R-squared:  0.6487 
F-statistic: 6.276 on 7 and 13 DF,  p-value: 0.002282

Resultados:

• Valor del coeficiente de determinación: 0.7717. Es un valor aceptable, el modelo únicamente explica el 77.17% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM2_Semanal <- forecast(fitSM2_Semanal, h=7, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM2_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Modelo

fitSM3_Semanal <- tslm(tsSM3_Semanal ~ trend + season) 
summary(fitSM3_Semanal)

Call:
tslm(formula = tsSM3_Semanal ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-50070 -20032  -1703  18539  49320 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   410909      20033  20.511 2.76e-11 ***
trend           5485       1171   4.684 0.000427 ***
season2         4110      25070   0.164 0.872309    
season3        -5833      25152  -0.232 0.820231    
season4       -45960      25288  -1.817 0.092260 .  
season5        16778      25477   0.659 0.521682    
season6        36018      25718   1.401 0.184776    
season7       -48786      26010  -1.876 0.083330 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 30670 on 13 degrees of freedom
Multiple R-squared:  0.7601,    Adjusted R-squared:  0.6309 
F-statistic: 5.884 on 7 and 13 DF,  p-value: 0.00305

Resultados:

• Valor del coeficiente de determinación: 0.7601. Es un valor que podría ser aceptable, el modelo únicamente explica el 76% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM3_Semanal <- forecast(fitSM3_Semanal, h=7, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM3_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Submetering 1

 tsSM1_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend            remainder         
 Min.   :-27570.39   Min.   :85332.42   Min.   :-17842.000  
 1st Qu.:-27284.26   1st Qu.:86567.92   1st Qu.: -2325.906  
 Median :-17054.28   Median :88045.73   Median :  -446.063  
 Mean   :     0.00   Mean   :88659.80   Mean   : -1989.231  
 3rd Qu.: 45539.85   3rd Qu.:90927.75   3rd Qu.:   621.535  
 Max.   : 59693.44   Max.   :91708.02   Max.   :  5531.746  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     72824         4360      2947         65904
   % 110.5          6.6       4.5         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.7882  0.9452  0.7781  0.9947  0.9999 

 Other components: List of 5
 $ win  : Named num [1:3] 211 11 7
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 22 2 1
 $ inner: int 1
 $ outer: int 15

tsSM1_Semanal %>%
  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
  autoplot()

Remainder no es un muelle. Observamos tendencia creciente y fuerte componente estacional

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

Componentes_SM1_SemanalSTL<- tsSM1_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM1_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM1_Semanal, col = 'gray')
lines(Componentes_SM1_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 2

 tsSM2_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE) %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend             remainder        
 Min.   :-39520.48   Min.   : 76016.18   Min.   :-20665.91  
 1st Qu.:-36421.48   1st Qu.: 78225.80   1st Qu.: -5262.30  
 Median : -4900.90   Median : 85975.97   Median :  -975.05  
 Mean   :     0.00   Mean   : 95499.46   Mean   :  5065.54  
 3rd Qu.: 29057.65   3rd Qu.:112604.73   3rd Qu.:  7002.00  
 Max.   : 47277.29   Max.   :133782.17   Max.   : 62291.18  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     65479        34379     12264         60848
   % 107.6         56.5      20.2         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8936  0.9452  0.8378  0.9882  0.9991 

 Other components: List of 5
 $ win  : Named num [1:3] 211 11 7
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 22 2 1
 $ inner: int 1
 $ outer: int 15

tsSM2_Semanal %>%
  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
  autoplot()

Componentes_SM2_SemanalSTL<- tsSM2_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM2_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM2_Semanal, col = 'gray')
lines(Componentes_SM2_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

Submetering 3

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia creciente

tsSM3_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend            remainder         
 Min.   :-46280.96   Min.   :444191.6   Min.   :-117714.08  
 1st Qu.:-44234.04   1st Qu.:449466.2   1st Qu.: -10712.44  
 Median : -5540.13   Median :456724.5   Median :  -6746.47  
 Mean   :     0.00   Mean   :470965.9   Mean   :  -5957.36  
 3rd Qu.: 40335.65   3rd Qu.:489708.8   3rd Qu.:  11362.30  
 Max.   : 41697.85   Max.   :534471.1   Max.   :  29107.15  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     84570        40243     22075         63070
   % 134.1         63.8      35.0         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.9035  0.9452  0.8755  0.9737  0.9937 

 Other components: List of 5
 $ win  : Named num [1:3] 211 11 7
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 22 2 1
 $ inner: int 1
 $ outer: int 15

Componentes_SM3_SemanalSTL<- tsSM3_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM3_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM3_Semanal, col = 'gray')
lines(Componentes_SM3_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Tendencia creciente del consumo energético

Evolución mensual del consumo (sin estacionalidad). Años 2007-2009

#tsSM1_Mensual_Adjusted <- tsSM1_Mensual - Componentes_SM1_Mensual$seasonal
#autoplot(tsSM1_Mensual_Adjusted)

tsSM1_Semanal_detrended <- tsSM1_Semanal - Componentes_SM1_SemanalSTL$time.series[,2]
plot.ts(tsSM1_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM1_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM1_SemanalSTL$time.series[,3])
qqline(Componentes_SM1_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM1_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM1_SemanalSTL$time.series[, 3]
W = 0.80542, p-value = 0.0007949

## Volvemos a descomponer la serie
plot(stl(tsSM1_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_Semanal <- HoltWinters(tsSM1_Semanal_detrended ,
                              #  tsSM1_Mensual_Adjusted, 
  beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Semanal_for <- forecast(tsSM1_HW_Semanal, h=7)
plot(tsSM1_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético semanal de la cocina \n Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Semanal_forC <- forecast(tsSM1_HW_Semanal, h=7, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción del consumo energético semanal para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM2_Mensual_Adjusted <- tsSM2_Mensual - Componentes_SM2_Mensual$seasonal
# autoplot(tsSM2_Mensual_Adjusted)

tsSM2_Semanal_detrended <- tsSM2_Semanal - Componentes_SM2_SemanalSTL$time.series[,2]
plot.ts(tsSM2_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM2_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM2_SemanalSTL$time.series[,3])
qqline(Componentes_SM2_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM2_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM2_SemanalSTL$time.series[, 3]
W = 0.70675, p-value = 3.345e-05

## Volvemos a descomponer la serie
plot(stl(tsSM2_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_Semanal <- HoltWinters(tsSM2_Semanal_detrended ,
                                # tsSM2_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
     tsSM2_HW_Semanal_for <- forecast(tsSM2_HW_Semanal, h=7)
plot(tsSM2_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción semanal del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Semanal_forC <- forecast(tsSM2_HW_Semanal, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción semanal para el año 2010", start(2010))

# ## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM3_Mensual_Adjusted <- tsSM3_Mensual - Componentes_SM3_Mensual$seasonal
# autoplot(tsSM3_Mensual_Adjusted)

tsSM3_Semanal_detrended <- tsSM3_Semanal - Componentes_SM3_SemanalSTL$time.series[,2]
plot.ts(tsSM3_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM3_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM3_SemanalSTL$time.series[,3])
qqline(Componentes_SM3_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM3_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM3_SemanalSTL$time.series[, 3]
W = 0.67378, p-value = 1.32e-05

## Volvemos a descomponer la serie
plot(stl(tsSM3_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_Semanal <- HoltWinters(tsSM3_Semanal_detrended ,
                              # tsSM3_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Semanal_for <- forecast(tsSM3_HW_Semanal, h=7)
plot(tsSM3_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético semanal del AC y termo eléctrico \n Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Semanal_forC <- forecast(tsSM3_HW_Semanal, h=7, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción semanal del consumo energético para el año 2010", start(2010))

RMSE_SM1_GranSemanalFor<-accuracy(tsSM1_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_1)[4]
RMSE_SM2_GranSemanalFor<-accuracy(tsSM2_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_2)[4]
RMSE_SM3_GranSemanalFor<-accuracy(tsSM3_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_3)[4]


ResEvolSemanal2_SinEstacionalidadConTendencia <-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranSemanalFor,
               RMSE_SM2_GranSemanalFor,
               RMSE_SM3_GranSemanalFor)
)
ResEvolSemanal2_SinEstacionalidadConTendencia %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Evolución mensual del consumo (con estacionalidad). Años 2007-2009

Vamos a ver una predicción del consumo enegético de cada submedidor con el método de HW. En este caso, no eliminaremos la estacionalidad, para ver así si se mantiene la tendencia de consumo.

tsSM1_Mensual_AdjustedII <- tsSM1_Mensual # - Componentes_SM1_Mensual$seasonal
autoplot(tsSM1_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM1_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_MensualII <- HoltWinters(tsSM1_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Mensual_forII <- forecast(tsSM1_HW_MensualII, h=12)
plot(tsSM1_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Mensual_forCII <- forecast(tsSM1_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM2_Mensual_AdjustedII <- tsSM2_Mensual ## - Componentes_SM2_Mensual$seasonal
autoplot(tsSM2_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM2_Mensual_AdjustedII)) 

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_MensualII <- HoltWinters(tsSM2_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_Mensual_forII <- forecast(tsSM2_HW_MensualII, h=12)
plot(tsSM2_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Mensual_forCII <- forecast(tsSM2_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM3_Mensual_AdjustedII <- tsSM3_Mensual # - Componentes_SM3_Mensual$seasonal
autoplot(tsSM3_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM3_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_MensualII <- HoltWinters(tsSM3_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Mensual_forII <- forecast(tsSM3_HW_MensualII, h=12)
plot(tsSM3_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Mensual_forCII <- forecast(tsSM3_HW_MensualII, h=12, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

RMSE_SM1_GranMensualForII<-accuracy(tsSM1_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_1)[4]
RMSE_SM2_GranMensualForII<-accuracy(tsSM2_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_2)[4]
RMSE_SM3_GranMensualForII<-accuracy(tsSM3_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_3)[4]


ResEvolMensual2SinTendenciaConEstacionalidad<-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranMensualForII,
               RMSE_SM2_GranMensualForII,
               RMSE_SM3_GranMensualForII)
   
    
)

ResEvolMensual2SinTendenciaConEstacionalidad %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Comparación de resultados. Predicción antes y despues de eliminar estacionalidad

Comparacion<-cbind.data.frame(
  ResEvolMensual2_SinEstacionalidadConTendencia,
  ResEvolMensual2SinTendenciaConEstacionalidad[,2])

colnames(Comparacion) <- c("Submedidor","RMSE Con Estacionalidad","RMSE Sin estacionalidad")
Comparacion %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Resultados más fiables y con menor error trás haber eliminado la estacionalidad de la serie. También así podemos ver la tendencia.

#p1<-autoplot(tsSM1_Mensual_Adjusted)
#p2<-autoplot(tsSM1_Mensual_AdjustedII)
#grid.arrange(
#             grobs = list(p1,p2),
#           nrow = 2,
#             top="Sin estacionalidad",bottom="Con estacionalidad"
#             )

Modelo

fitSM1_Semanal <- tslm(tsSM1_Semanal ~ trend + season) 
summary(fitSM1_Semanal)

Call:
tslm(formula = tsSM1_Semanal ~ trend + season)

Residuals:
     Min       1Q   Median       3Q      Max 
-10942.7  -1970.9   -182.3   2205.3  12203.0 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 10835.442   2693.122   4.023 0.000108 ***
trend          -8.872      8.426  -1.053 0.294745    
season2      6872.206   3753.914   1.831 0.069985 .  
season3      6643.745   3753.942   1.770 0.079662 .  
season4      4273.617   3753.989   1.138 0.257536    
season5      3290.489   3754.056   0.877 0.382751    
season6      4326.029   3754.141   1.152 0.251800    
season7       791.901   3754.245   0.211 0.833347    
season8      -313.227   3754.368  -0.083 0.933669    
season9     -3984.688   3754.509  -1.061 0.290985    
season10     6229.851   3754.670   1.659 0.100054    
season11     2045.390   3754.850   0.545 0.587092    
season12     5559.929   3755.048   1.481 0.141693    
season13     3915.802   3755.266   1.043 0.299459    
season14     1586.341   3755.502   0.422 0.673594    
season15     3100.547   3755.757   0.826 0.410934    
 [ reached getOption("max.print") -- omitted 38 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4598 on 105 degrees of freedom
Multiple R-squared:  0.4836,    Adjusted R-squared:  0.223 
F-statistic: 1.855 on 53 and 105 DF,  p-value: 0.003663

Resultados:

• Valor del coeficiente de determinación: 0.4836. Es un valor bajo, el modelo únicamente explica el 48% de la variabilidad total.

• Hay un coeficiente significativamente no nulo.

Predicción IC

Predicción para los próximos 20 días

forecastfitSM1_Semanal <- forecast(fitSM1_Semanal, h=53, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM1_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético semanal total del año 2010 \n  a través de un modelo de regresión")

Modelo

fitSM2_Semanal <- tslm(tsSM2_Semanal ~ trend + season) 
summary(fitSM2_Semanal)

Call:
tslm(formula = tsSM2_Semanal ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-12710  -2492    262   2070   9448 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)  19974.287   2851.564   7.005 2.45e-10 ***
trend          -46.079      8.921  -5.165 1.15e-06 ***
season2      -4864.254   3974.765  -1.224 0.223773    
season3        714.159   3974.795   0.180 0.857756    
season4       -986.429   3974.845  -0.248 0.804489    
season5      -1733.682   3974.915  -0.436 0.663619    
season6       -145.270   3975.005  -0.037 0.970917    
season7      -4670.524   3975.115  -1.175 0.242676    
season8      -1727.111   3975.245  -0.434 0.664841    
season9      -5223.698   3975.395  -1.314 0.191707    
season10      3920.048   3975.565   0.986 0.326383    
season11      1273.127   3975.756   0.320 0.749436    
season12     -1349.793   3975.966  -0.339 0.734920    
season13      -436.381   3976.196  -0.110 0.912818    
season14     -2924.635   3976.446  -0.735 0.463681    
season15     -2526.222   3976.716  -0.635 0.526645    
 [ reached getOption("max.print") -- omitted 38 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 4868 on 105 degrees of freedom
Multiple R-squared:  0.5326,    Adjusted R-squared:  0.2966 
F-statistic: 2.257 on 53 and 105 DF,  p-value: 0.0001984

Resultados:

• Valor del coeficiente de determinación: 0.5326. Es un valor aceptable, el modelo únicamente explica el 53% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM2_Semanal <- forecast(fitSM2_Semanal, h=53, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM2_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Modelo

fitSM3_Semanal <- tslm(tsSM3_Semanal ~ trend + season) 
summary(fitSM3_Semanal)

Call:
tslm(formula = tsSM3_Semanal ~ trend + season)

Residuals:
   Min     1Q Median     3Q    Max 
-35933  -5644    -17   6305  26151 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  56824.99    7005.03   8.112 9.94e-13 ***
trend           95.69      21.92   4.366 2.97e-05 ***
season2      13899.65    9764.23   1.424  0.15755    
season3      18434.30    9764.31   1.888  0.06180 .  
season4       9370.94    9764.43   0.960  0.33941    
season5      20968.59    9764.60   2.147  0.03406 *  
season6      14556.57    9764.82   1.491  0.13903    
season7       8915.22    9765.09   0.913  0.36335    
season8      -2728.46    9765.41  -0.279  0.78049    
season9     -25352.82    9765.78  -2.596  0.01078 *  
season10      2945.17    9766.20   0.302  0.76358    
season11      8872.81    9766.67   0.908  0.36571    
season12      8371.46    9767.18   0.857  0.39334    
season13     13322.78    9767.75   1.364  0.17550    
season14      9586.09    9768.36   0.981  0.32868    
season15     -1125.59    9769.03  -0.115  0.90849    
 [ reached getOption("max.print") -- omitted 38 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 11960 on 105 degrees of freedom
Multiple R-squared:  0.7153,    Adjusted R-squared:  0.5715 
F-statistic: 4.977 on 53 and 105 DF,  p-value: 1.341e-12

Resultados:

• Valor del coeficiente de determinación: 0.7153. Es un valor que podría ser aceptable, el modelo únicamente explica el 71% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM3_Semanal <- forecast(fitSM3_Semanal, h=53, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM3_Semanal, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Submetering 1

 tsSM1_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal              trend             remainder         
 Min.   :-11380.560   Min.   : 9373.389   Min.   :-17306.164  
 1st Qu.: -1525.527   1st Qu.:10411.060   1st Qu.: -1565.107  
 Median :   605.971   Median :11525.838   Median :    50.446  
 Mean   :     0.000   Mean   :11285.767   Mean   :   161.290  
 3rd Qu.:  2468.518   3rd Qu.:12047.600   3rd Qu.:  1583.123  
 Max.   : 12092.727   Max.   :12804.151   Max.   : 18204.879  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     3994         1637      3148          6386 
   %  62.5         25.6      49.3         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8303  0.9452  0.8244  0.9824  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 1591 81 53
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 160 9 6
 $ inner: int 1
 $ outer: int 15

#tsSM1_Semanal %>%
#  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
#  autoplot()

Remainder no es un muelle. Observamos tendencia creciente y fuerte componente estacional

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

Componentes_SM1_SemanalSTL<- tsSM1_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM1_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM1_Semanal, col = 'gray')
lines(Componentes_SM1_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 2

 tsSM2_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE) %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal              trend            remainder         
 Min.   :-11082.099   Min.   :10972.07   Min.   :-17744.999  
 1st Qu.: -2312.887   1st Qu.:11343.85   1st Qu.: -1509.064  
 Median :   306.932   Median :12813.79   Median :   128.053  
 Mean   :     0.000   Mean   :13181.15   Mean   :   101.019  
 3rd Qu.:  2737.120   3rd Qu.:15219.46   3rd Qu.:  1570.218  
 Max.   :  9565.068   Max.   :15803.51   Max.   : 14301.230  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     5050         3876      3079          7916 
   %  63.8         49.0      38.9         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.7600  0.9452  0.7848  0.9927  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 1591 81 53
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 160 9 6
 $ inner: int 1
 $ outer: int 15

#tsSM2_Semanal %>%
#  stl(t.window=13, s.window="periodic", robust=TRUE) %>%
#  autoplot()

Componentes_SM2_SemanalSTL<- tsSM2_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM2_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM2_Semanal, col = 'gray')
lines(Componentes_SM2_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

Submetering 3

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia creciente

tsSM3_Semanal %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend            remainder        
 Min.   :-47211.03   Min.   :57753.65   Min.   :-56145.14  
 1st Qu.:-10494.71   1st Qu.:58324.34   1st Qu.: -4929.51  
 Median :  4528.45   Median :63141.11   Median :  -321.52  
 Mean   :     0.00   Mean   :62123.83   Mean   :  -707.60  
 3rd Qu.: 11265.72   3rd Qu.:63648.44   3rd Qu.:  3997.01  
 Max.   : 26621.30   Max.   :68148.68   Max.   : 39973.53  
 IQR:
     STL.seasonal STL.trend STL.remainder data 
     21760         5324      8927         23820
   %  91.4         22.4      37.5         100.0

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8598  0.9452  0.8486  0.9898  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 1591 81 53
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 160 9 6
 $ inner: int 1
 $ outer: int 15

Componentes_SM3_SemanalSTL<- tsSM3_Semanal %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM3_SemanalSTL, ts.colour = 'blue')

plot.ts(tsSM3_Semanal, col = 'gray')
lines(Componentes_SM3_SemanalSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Tendencia creciente del consumo energético

Evolución mensual del consumo (sin estacionalidad). Años 2007-2009

#tsSM1_Mensual_Adjusted <- tsSM1_Mensual - Componentes_SM1_Mensual$seasonal
#autoplot(tsSM1_Mensual_Adjusted)

tsSM1_Semanal_detrended <- tsSM1_Semanal - Componentes_SM1_SemanalSTL$time.series[,2]
plot.ts(tsSM1_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM1_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM1_SemanalSTL$time.series[,3])
qqline(Componentes_SM1_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM1_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM1_SemanalSTL$time.series[, 3]
W = 0.83235, p-value = 3.231e-12

## Volvemos a descomponer la serie
plot(stl(tsSM1_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_Semanal <- HoltWinters(tsSM1_Semanal_detrended ,
                              #  tsSM1_Mensual_Adjusted, 
  beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Semanal_for <- forecast(tsSM1_HW_Semanal, h=53)
plot(tsSM1_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético semanal de la cocina \n Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Semanal_forC <- forecast(tsSM1_HW_Semanal, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción del consumo energético semanal para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM2_Mensual_Adjusted <- tsSM2_Mensual - Componentes_SM2_Mensual$seasonal
# autoplot(tsSM2_Mensual_Adjusted)

tsSM2_Semanal_detrended <- tsSM2_Semanal - Componentes_SM2_SemanalSTL$time.series[,2]
plot.ts(tsSM2_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM2_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM2_SemanalSTL$time.series[,3])
qqline(Componentes_SM2_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM2_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM2_SemanalSTL$time.series[, 3]
W = 0.88071, p-value = 5.464e-10

## Volvemos a descomponer la serie
plot(stl(tsSM2_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_Semanal <- HoltWinters(tsSM2_Semanal_detrended ,
                                # tsSM2_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
     tsSM2_HW_Semanal_for <- forecast(tsSM2_HW_Semanal, h=53)
plot(tsSM2_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción semanal del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Semanal_forC <- forecast(tsSM2_HW_Semanal, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción semanal para el año 2010", start(2010))

# ## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
# tsSM3_Mensual_Adjusted <- tsSM3_Mensual - Componentes_SM3_Mensual$seasonal
# autoplot(tsSM3_Mensual_Adjusted)

tsSM3_Semanal_detrended <- tsSM3_Semanal - Componentes_SM3_SemanalSTL$time.series[,2]
plot.ts(tsSM3_Semanal_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
plot(Componentes_SM3_SemanalSTL$time.series[,3], 
     col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM3_SemanalSTL$time.series[,3])
qqline(Componentes_SM3_SemanalSTL$time.series[,3])

shapiro.test(Componentes_SM3_SemanalSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM3_SemanalSTL$time.series[, 3]
W = 0.83537, p-value = 4.311e-12

## Volvemos a descomponer la serie
plot(stl(tsSM3_Semanal_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_Semanal <- HoltWinters(tsSM3_Semanal_detrended ,
                              # tsSM3_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_Semanal,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Semanal_for <- forecast(tsSM3_HW_Semanal, h=53)
plot(tsSM3_HW_Semanal_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético semanal del AC y termo eléctrico \n Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Semanal_forC <- forecast(tsSM3_HW_Semanal, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Semanal_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción semanal del consumo energético para el año 2010", start(2010))

RMSE_SM1_GranSemanalFor<-accuracy(tsSM1_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_1)[4]
RMSE_SM2_GranSemanalFor<-accuracy(tsSM2_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_2)[4]
RMSE_SM3_GranSemanalFor<-accuracy(tsSM3_HW_Semanal_forC   ,DatosRealesSeriesSemanales$Sub_metering_3)[4]


ResEvolSemanal2_SinEstacionalidadConTendencia <-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranSemanalFor,
               RMSE_SM2_GranSemanalFor,
               RMSE_SM3_GranSemanalFor)
)
ResEvolSemanal2_SinEstacionalidadConTendencia %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Evolución mensual del consumo (con estacionalidad). Años 2007-2009

Vamos a ver una predicción del consumo enegético de cada submedidor con el método de HW. En este caso, no eliminaremos la estacionalidad, para ver así si se mantiene la tendencia de consumo.

tsSM1_Mensual_AdjustedII <- tsSM1_Mensual # - Componentes_SM1_Mensual$seasonal
autoplot(tsSM1_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM1_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_MensualII <- HoltWinters(tsSM1_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_Mensual_forII <- forecast(tsSM1_HW_MensualII, h=12)
plot(tsSM1_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_Mensual_forCII <- forecast(tsSM1_HW_MensualII, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM2_Mensual_AdjustedII <- tsSM2_Mensual ## - Componentes_SM2_Mensual$seasonal
autoplot(tsSM2_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM2_Mensual_AdjustedII)) 

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_MensualII <- HoltWinters(tsSM2_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_Mensual_forII <- forecast(tsSM2_HW_MensualII, h=53)
plot(tsSM2_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_Mensual_forCII <- forecast(tsSM2_HW_MensualII, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM3_Mensual_AdjustedII <- tsSM3_Mensual # - Componentes_SM3_Mensual$seasonal
autoplot(tsSM3_Mensual_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM3_Mensual_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_MensualII <- HoltWinters(tsSM3_Mensual_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_MensualII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_Mensual_forII <- forecast(tsSM3_HW_MensualII, h=53)
plot(tsSM3_HW_Mensual_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_Mensual_forCII <- forecast(tsSM3_HW_MensualII, h=53, level=c(10,25))
## Plot only the forecasted area
plot(tsSM3_HW_Mensual_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

RMSE_SM1_GranMensualForII<-accuracy(tsSM1_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_1)[4]
RMSE_SM2_GranMensualForII<-accuracy(tsSM2_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_2)[4]
RMSE_SM3_GranMensualForII<-accuracy(tsSM3_HW_Mensual_forCII   ,DatosRealesSeriesMensuales$Sub_metering_3)[4]


ResEvolMensual2SinTendenciaConEstacionalidad<-
cbind.data.frame(
      Submedidor=c("Cocina","Lavadero","AC y TermoE"),
      RMSE = c(RMSE_SM1_GranMensualForII,
               RMSE_SM2_GranMensualForII,
               RMSE_SM3_GranMensualForII)
   
    
)

ResEvolMensual2SinTendenciaConEstacionalidad %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Comparación de resultados. Predicción antes y despues de eliminar estacionalidad

Comparacion<-cbind.data.frame(
  ResEvolMensual2_SinEstacionalidadConTendencia,
  ResEvolMensual2SinTendenciaConEstacionalidad[,2])

colnames(Comparacion) <- c("Submedidor","RMSE Con Estacionalidad","RMSE Sin estacionalidad")
Comparacion %>% kable(booktabs=TRUE) %>% kable_styling(latex_options = "striped")

Resultados más fiables y con menor error trás haber eliminado la estacionalidad de la serie. También así podemos ver la tendencia.

#p1<-autoplot(tsSM1_Mensual_Adjusted)
#p2<-autoplot(tsSM1_Mensual_AdjustedII)
#grid.arrange(
#             grobs = list(p1,p2),
#           nrow = 2,
#             top="Sin estacionalidad",bottom="Con estacionalidad"
#             )

Modelo

fitSM1_DiariaEnero <- tslm(tsSM1_DiariaEnero ~ trend + season) 
summary(fitSM1_DiariaEnero)

Call:
tslm(formula = tsSM1_DiariaEnero ~ trend + season)

Residuals:
    Min      1Q  Median      3Q     Max 
-3996.9 -1161.0  -265.1   822.5  6806.5 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)  
(Intercept)  771.999   1295.981   0.596   0.5536  
trend          9.479      8.968   1.057   0.2947  
season2     -279.146   1787.319  -0.156   0.8764  
season3     1429.042   1787.386   0.800   0.4271  
season4      145.896   1787.499   0.082   0.9352  
season5      501.083   1787.656   0.280   0.7802  
season6      -22.063   1787.859  -0.012   0.9902  
season7      855.125   1788.106   0.478   0.6342  
season8      514.979   1788.399   0.288   0.7744  
season9      399.500   1788.736   0.223   0.8240  
season10     474.687   1789.118   0.265   0.7917  
season11     865.541   1789.545   0.484   0.6304  
season12    1036.062   1790.017   0.579   0.5649  
season13    3653.916   1790.533   2.041   0.0456 *
season14    2432.437   1791.095   1.358   0.1794  
season15     330.291   1791.701   0.184   0.8544  
 [ reached getOption("max.print") -- omitted 16 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2189 on 61 degrees of freedom
Multiple R-squared:  0.2234,    Adjusted R-squared:  -0.1713 
F-statistic: 0.566 on 31 and 61 DF,  p-value: 0.957

Resultados:

• Valor del coeficiente de determinación: 0.2234. Es un valor muy bajo, el modelo únicamente explica el 22% de la variabilidad total.

• Hay un coeficiente significativamente no nulo.

Predicción IC

Predicción para el mes de Enero de 2010

forecastfitSM1_DiariaEnero <- forecast(fitSM1_DiariaEnero, h=31, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM1_DiariaEnero, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético de Enero del año 2010 \n  a través de un modelo de regresión")

Modelo

fitSM2_DiariaEnero <- tslm(tsSM2_DiariaEnero ~ trend + season) 
summary(fitSM2_DiariaEnero)

Call:
tslm(formula = tsSM2_DiariaEnero ~ trend + season)

Residuals:
    Min      1Q  Median      3Q     Max 
-4045.7 -1317.7  -282.9  1418.7  4220.7 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)  
(Intercept)  2217.675   1483.490   1.495   0.1401  
trend          -3.261     10.265  -0.318   0.7518  
season2       696.594   2045.917   0.340   0.7347  
season3     -1074.479   2045.995  -0.525   0.6014  
season4      3980.449   2046.123   1.945   0.0563 .
season5     -1766.291   2046.304  -0.863   0.3914  
season6     -1070.697   2046.535  -0.523   0.6027  
season7      1986.897   2046.819   0.971   0.3355  
season8      -674.175   2047.153  -0.329   0.7430  
season9      -885.581   2047.539  -0.433   0.6669  
season10      900.346   2047.977   0.440   0.6618  
season11      927.940   2048.465   0.453   0.6522  
season12     -943.799   2049.005  -0.461   0.6467  
season13     -847.205   2049.597  -0.413   0.6808  
season14     -450.611   2050.239  -0.220   0.8268  
season15      594.649   2050.933   0.290   0.7728  
 [ reached getOption("max.print") -- omitted 16 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2506 on 61 degrees of freedom
Multiple R-squared:  0.2859,    Adjusted R-squared:  -0.07693 
F-statistic: 0.788 on 31 and 61 DF,  p-value: 0.7627

Resultados:

• Valor del coeficiente de determinación: 0.2859. Es un valor bajo, el modelo únicamente explica el 29% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM2_DiariaEnero <- forecast(fitSM2_DiariaEnero, h=31, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM2_DiariaEnero, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del año 2010 \n  a través de un modelo de regresión")

Modelo

fitSM3_DiariaEnero <- tslm(tsSM3_DiariaEnero~ trend + season) 
summary(fitSM3_DiariaEnero)

Call:
tslm(formula = tsSM3_DiariaEnero ~ trend + season)

Residuals:
    Min      1Q  Median      3Q     Max 
-5223.2 -1643.5  -128.2  1968.5  6090.7 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 7359.3708  1811.1381   4.063  0.00014 ***
trend         -0.8762    12.5325  -0.070  0.94449    
season2     1941.2095  2497.7856   0.777  0.44006    
season3     -337.9143  2497.8799  -0.135  0.89284    
season4     1701.2952  2498.0371   0.681  0.49842    
season5     2241.8380  2498.2572   0.897  0.37305    
season6     1754.3809  2498.5401   0.702  0.48525    
season7     3264.9237  2498.8858   1.307  0.19627    
season8     5867.7999  2499.2943   2.348  0.02215 *  
season9     2021.6760  2499.7656   0.809  0.42180    
season10    3640.5522  2500.2996   1.456  0.15051    
season11    1613.7617  2500.8963   0.645  0.52117    
season12    2887.3045  2501.5556   1.154  0.25292    
season13    4624.5140  2502.2776   1.848  0.06943 .  
season14    3936.0569  2503.0621   1.572  0.12101    
season15    5955.9331  2503.9090   2.379  0.02052 *  
 [ reached getOption("max.print") -- omitted 16 rows ]
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3059 on 61 degrees of freedom
Multiple R-squared:  0.3269,    Adjusted R-squared:  -0.01514 
F-statistic: 0.9557 on 31 and 61 DF,  p-value: 0.5438

Resultados:

• Valor del coeficiente de determinación:0.3269. Es un valor que es muy bajo, el modelo únicamente explica el 32% de la variabilidad total.

• Hay tres coeficientes significativamente no nulos.

Predicción IC

Predicción para los próximos 12 meses (año 2010)

forecastfitSM3_DiariaEnero <- forecast(fitSM3_DiariaEnero, h=31, level=c(80,90))
## h=12 para predecir el año 2010
## Plot sub-meter 3 forecast, limit y and add labels
plot(forecastfitSM3_DiariaEnero, ylab= "Watt-Hours", xlab="Fecha",
     main="Predicción del consumo energético del mes de Enero del año 2010 \n  a través de un modelo de regresión")

Se prevee una bajada del consumo

Submetering 1

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia decreciente.

 tsSM1_DiariaEnero %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend             remainder        
 Min.   :-1324.873   Min.   : 871.8631   Min.   :-6121.092  
 1st Qu.: -761.466   1st Qu.:1309.3572   1st Qu.: -491.953  
 Median : -158.424   Median :1491.6212   Median :  220.447  
 Mean   :    0.000   Mean   :1522.7948   Mean   :  550.818  
 3rd Qu.:  337.611   3rd Qu.:1775.0069   3rd Qu.:  732.488  
 Max.   : 5397.352   Max.   :2077.1908   Max.   : 9934.040  
 IQR:
     STL.seasonal STL.trend STL.remainder data  
     1099.1        465.6    1224.4        1692.0
   %  65.0         27.5      72.4         100.0 

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8124  0.9452  0.7973  0.9880  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 931 47 31
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 94 5 4
 $ inner: int 1
 $ outer: int 15

Remainder no es un muelle. Observamos tendencia creciente y fuerte componente estacional

Componentes_SM1_DiariaEneroSTL<-
  tsSM1_DiariaEnero %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM1_DiariaEneroSTL, ts.colour = 'blue')

plot.ts(tsSM1_Mensual, col = 'gray')
lines(Componentes_SM1_DiariaEneroSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 2

 tsSM2_DiariaEnero %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE) %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend            remainder        
 Min.   :-1788.921   Min.   :1491.391   Min.   :-5599.893  
 1st Qu.:-1294.941   1st Qu.:1617.887   1st Qu.: -685.829  
 Median : -782.743   Median :1866.286   Median : -122.576  
 Mean   :    0.000   Mean   :1939.768   Mean   :  373.264  
 3rd Qu.:  718.901   3rd Qu.:2339.704   3rd Qu.: 1313.902  
 Max.   : 4203.925   Max.   :2437.168   Max.   : 6567.217  
 IQR:
     STL.seasonal STL.trend STL.remainder data  
     2013.8        721.8    1999.7        3494.0
   %  57.6         20.7      57.2         100.0 

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8067  0.9452  0.7789  0.9886  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 931 47 31
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 94 5 4
 $ inner: int 1
 $ outer: int 15

Componentes_SM2_DiariaEneroSTL<- tsSM2_DiariaEnero %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM2_DiariaEneroSTL, ts.colour = 'blue')

plot.ts(tsSM2_DiariaEnero, col = 'gray')
lines(Componentes_SM2_DiariaEneroSTL$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Submetering 3

## Descomponer la serie
Componentes_SM3_DiariaEnero <- decompose(tsSM3_DiariaEnero)
## Plot
plot(Componentes_SM3_DiariaEnero)

## Check summary statistics for decomposed sub-meter 2
summary(Componentes_SM3_DiariaEnero)

         Length Class  Mode     
x        93     ts     numeric  
seasonal 93     ts     numeric  
trend    93     ts     numeric  
random   93     ts     numeric  
figure   31     -none- numeric  
type      1     -none- character

Se aprecia componente estacional: oscilacione que se repiten cada año. La componente tendencia tambien se aprecia. Tendencia creciente

tsSM3_DiariaEnero %>%
  stl( #t.window=13,
      s.window="periodic", 
      robust=TRUE)  %>% summary()

 Call:
 stl(x = ., s.window = "periodic", robust = TRUE)

 Time.series components:
    seasonal             trend             remainder        
 Min.   :-5596.817   Min.   : 9448.476   Min.   :-8001.466  
 1st Qu.:-1429.752   1st Qu.:10055.012   1st Qu.:-1298.672  
 Median :  135.662   Median :10467.966   Median : -172.626  
 Mean   :    0.000   Mean   :10408.787   Mean   :   31.439  
 3rd Qu.: 1423.036   3rd Qu.:10742.199   3rd Qu.: 1502.564  
 Max.   : 3704.229   Max.   :11008.213   Max.   : 8713.785  
 IQR:
     STL.seasonal STL.trend STL.remainder data  
     2852.8        687.2    2801.2        3550.0
   %  80.4         19.4      78.9         100.0 

 Weights:
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.8394  0.9452  0.8553  0.9931  1.0000 

 Other components: List of 5
 $ win  : Named num [1:3] 931 47 31
 $ deg  : Named int [1:3] 0 1 1
 $ jump : Named num [1:3] 94 5 4
 $ inner: int 1
 $ outer: int 15

Componentes_SM3_DiariaEneroSTL<- tsSM3_DiariaEnero %>%stl( #t.window=13,
  s.window="periodic",  robust=TRUE)
autoplot(Componentes_SM3_DiariaEneroSTL, ts.colour = 'blue')

plot.ts(tsSM3_DiariaEnero, col = 'gray')
lines(Componentes_SM3_DiariaEnero$time.series[,2], 
      col = "red", lwd = 1, lty = 2)

Evolución mensual del consumo (sin estacionalidad). Años 2007-2009

tsSM1_DiariaEnero_detrended <- tsSM1_DiariaEnero - Componentes_SM1_DiariaEneroSTL$time.series[,2]
plot.ts(tsSM1_DiariaEnero_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
  plot(Componentes_SM1_DiariaEneroSTL$time.series[,3], col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM1_DiariaEneroSTL$time.series[,3])
qqline(Componentes_SM1_DiariaEneroSTL$time.series[,3])

shapiro.test(Componentes_SM1_DiariaEneroSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM1_DiariaEneroSTL$time.series[, 3]
W = 0.78722, p-value = 2.838e-10

## Volvemos a descomponer la serie
plot(stl(tsSM1_DiariaEnero_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_DiariaEnero<- HoltWinters(tsSM1_DiariaEnero_detrended ,
                              #  tsSM1_Mensual_Adjusted, 
  beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_DiariaEnero,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_DiariaEnero_for <- forecast(tsSM1_HW_DiariaEnero, h=31)
plot(tsSM1_HW_DiariaEnero_for , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_DiariaEnero_forC <- forecast(tsSM1_HW_DiariaEnero, h=31, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_DiariaEnero_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

tsSM2_DiariaEnero_detrended <- tsSM2_DiariaEnero - Componentes_SM2_DiariaEneroSTL$time.series[,2]
plot.ts(tsSM2_DiariaEnero_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
  plot(Componentes_SM2_DiariaEneroSTL$time.series[,3], col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM2_DiariaEneroSTL$time.series[,3])
qqline(Componentes_SM2_DiariaEneroSTL$time.series[,3])

shapiro.test(Componentes_SM2_DiariaEneroSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM2_DiariaEneroSTL$time.series[, 3]
W = 0.88372, p-value = 5.75e-07

## Volvemos a descomponer la serie
plot(stl(tsSM2_DiariaEnero_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_DiariaEnero <- HoltWinters(tsSM2_DiariaEnero_detrended ,
                                # tsSM2_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_DiariaEnero,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_DiariaEnero_for <- forecast(tsSM2_HW_DiariaEnero, h=31)
plot(tsSM2_HW_DiariaEnero_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_DiariaEnero_forC <- forecast(tsSM2_HW_DiariaEnero, h=31, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_DiariaEnero_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

tsSM3_DiariaEnero_detrended <- tsSM3_DiariaEnero - Componentes_SM3_DiariaEneroSTL$time.series[,2]
plot.ts(tsSM3_DiariaEnero_detrended, 
        ylab = "Energía (Vatios-hora)", 
        main = 'Seasonal variability of energy comsumption',
        cex.main = 0.85)

par(mfrow = c(1,2))
  plot(Componentes_SM3_DiariaEneroSTL$time.series[,3], col = "blue", main = 'Remainder', ylab = "")
qqnorm(Componentes_SM3_DiariaEneroSTL$time.series[,3])
qqline(Componentes_SM3_DiariaEneroSTL$time.series[,3])

shapiro.test(Componentes_SM3_DiariaEneroSTL$time.series[,3])

    Shapiro-Wilk normality test

data:  Componentes_SM3_DiariaEneroSTL$time.series[, 3]
W = 0.95542, p-value = 0.003018

## Volvemos a descomponer la serie
plot(stl(tsSM3_DiariaEnero_detrended ,s.window = "periodic",robust = TRUE))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_DiariaEnero <- HoltWinters(tsSM3_DiariaEnero_detrended ,
                              # tsSM3_Mensual_Adjusted,
                                beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_DiariaEnero,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM3_HW_DiariaEnero_for <- forecast(tsSM3_HW_DiariaEnero, h=31)
plot(tsSM3_HW_DiariaEnero_for , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM3_HW_DiariaEnero_forC <- forecast( object = tsSM3_HW_DiariaEnero_for, h=31)
## Plot only the forecasted area
plot(tsSM3_HW_DiariaEnero_forC, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

# RMSE_SM1_GranDiariaEneroFor<-accuracy(tsSM1_HW_DiariaEnero_forC   # ,seriesDiariasEnero$Sub_metering_1)[4]
# RMSE_SM2_GranDiariaEneroFor<-accuracy(tsSM2_HW_DiariaEnero_forC   # ,seriesDiariasEnero$Sub_metering_2)[4]
# RMSE_SM3_GranDiariaEneroFor<-accuracy(tsSM3_HW_DiariaEnero_forC   # ,seriesDiariasEnero$Sub_metering_3)[4]
# 
# 
# ResEvolDiariaEnero2_SinEstacionalidadConTendencia <-
# cbind.data.frame(
#       Submedidor=c("Cocina","Lavadero","AC y TermoE"),
#       RMSE = c(RMSE_SM1_GranDiariaEneroFor,
#                RMSE_SM2_GranDiariaEneroFor,
#                RMSE_SM3_GranDiariaEneroFor)
# )
# ResEvolDiariaEnero2_SinEstacionalidadConTendencia %>% kable(booktabs=TRUE) %>% # kable_styling(latex_options = "striped")

Evolución mensual del consumo (con estacionalidad). Años 2007-2009

Vamos a ver una predicción del consumo enegético de cada submedidor con el método de HW. En este caso, no eliminaremos la estacionalidad, para ver así si se mantiene la tendencia de consumo.

tsSM1_DiariaEnero_AdjustedII <- tsSM1_DiariaEnero # - Componentes_SM1_Mensual$seasonal
autoplot(tsSM1_DiariaEnero_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM1_DiariaEnero_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM1_HW_DiariaEneroII <- HoltWinters(tsSM1_DiariaEnero_AdjustedII, beta=TRUE, gamma=TRUE)
# optimiza
plot(tsSM1_HW_DiariaEneroII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM1_HW_DiariaEnero_forII <- forecast(tsSM1_HW_DiariaEneroII, h=31)
plot(tsSM1_HW_DiariaEnero_forII , ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico",
     main="Predicción del consumo energético de la cocina \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM1_HW_DiariaEnero_forCII <- forecast(tsSM1_HW_DiariaEneroII, h=31, level=c(10,25))
## Plot only the forecasted area
plot(tsSM1_HW_DiariaEnero_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM2_DiariaEnero_AdjustedII <- tsSM2_DiariaEnero ## - Componentes_SM2_Mensual$seasonal
autoplot(tsSM2_DiariaEnero_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM2_DiariaEnero_AdjustedII)) 

## Holt Winters Exponential Smoothing & Plot
tsSM2_HW_DiariaEneroII <- HoltWinters(tsSM2_DiariaEnero_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM2_HW_DiariaEneroII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
tsSM2_HW_DiariaEnero_forII <- forecast(tsSM2_HW_DiariaEneroII, h=31)
plot(tsSM2_HW_DiariaEnero_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del lavadero \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
tsSM2_HW_DiariaEnero_forCII <- forecast(tsSM2_HW_DiariaEneroII, h=31, level=c(10,25))
## Plot only the forecasted area
plot(tsSM2_HW_DiariaEnero_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

## Seasonal adjusting sub-meter 3 by subtracting the seasonal component & plot
tsSM3_DiariaEnero_AdjustedII <- tsSM3_DiariaEnero # - Componentes_SM3_Mensual$seasonal
autoplot(tsSM3_DiariaEnero_AdjustedII)

## Volvemos a descomponer la serie
plot(decompose(tsSM3_DiariaEnero_AdjustedII))

## Holt Winters Exponential Smoothing & Plot
tsSM3_HW_DiariaEneroII <- HoltWinters(tsSM3_DiariaEnero_AdjustedII, beta=TRUE, gamma=TRUE)
plot(tsSM3_HW_DiariaEneroII,
     xlab="Fecha",
     ylab="Valores observados y ajustados",
     main="Suavizado exponencial de Holt-Winters")

## HoltWinters forecast & plot
     tsSM3_HW_DiariaEnero_forII <- forecast(tsSM3_HW_DiariaEneroII, h=31)
plot(tsSM3_HW_DiariaEnero_forII , ylab= "Vatios-Hora", xlab="Fecha - Lavadero",
     main="Predicción del consumo energético del AC y termo eléctrico \ Año 2010" ) 

## Forecast HoltWinters with diminished confidence levels
     tsSM3_HW_DiariaEnero_forCII <- forecast(tsSM3_HW_DiariaEneroII, h=31, level=c(10,25))
plot(tsSM3_HW_DiariaEnero_forCII, ylab= "Vatios-Hora", xlab="Fecha - AireAcondicionado y termo eléctrico", main="Predicción para el año 2010", start(2010))

# RMSE_SM1_GranDiariaEneroForII<-accuracy(tsSM1_HW_DiariaEnero_forCII   # ,seriesDiariasEnero$Sub_metering_1)[4]
# RMSE_SM2_GranDiariaEneroForII<-accuracy(tsSM2_HW_DiariaEnero_forCII   # ,seriesDiariasEnero$Sub_metering_2)[4]
# RMSE_SM3_GranDiariaEneroForII<-accuracy(tsSM3_HW_DiariaEnero_forCII   # ,seriesDiariasEnero$Sub_metering_3)[4]
# 
# 
# ResEvolDiariaEnero2SinTendenciaConEstacionalidad<-
# cbind.data.frame(
#       Submedidor=c("Cocina","Lavadero","AC y TermoE"),
#       RMSE = c(RMSE_SM1_GranDiariaEneroForII,
#                RMSE_SM2_GranDiariaEneroForII,
#                RMSE_SM3_GranDiariaEneroForII)
#    
#     
# )
# 
# ResEvolDiariaEnero2SinTendenciaConEstacionalidad %>% kable(booktabs=TRUE) %>% # kable_styling(latex_options = "striped")